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ABSTRACT 

We study the fully non-adiabatic tidal response of a uniformly rotating unevolved 
20Mq star to the dominant I = m = 2 component of the companion’s perturbing 
potential. This is done numerically with a 2D implicit finite difference scheme. We 
assume the star is rotating slowly with angular speed Dg flc, so that the centrifugal 
force can be neglected, but the Coriolis force is taken fully into account. We study the 
low frequency ‘inertial’ regime |ct| < 2ns, where a is the forcing frequency in the frame 
rotating with the stellar spin rate fig. In this frequency range inertial modes are excited 
in the convective core which can interact with rotationally modified g- or r-modes 
in the radiative envelope and cause significant strengthening of the tidal interaction. 
Resonant interaction with quasi-toroidal (r-)modes in slightly super-synchronous stars 
causes efficient spin down towards corotation. We determine timescales for tidal spin- 
up and spin-down in the inertial frequency regime for stars spinning with Dg = O.lDc 
and 0.2 Dc. 

Key words: Hydrodynamics- Stars: binaries- Stars: rotation- Stars: oscillation - 
Stars: tides 


1 INTRODUCTION 

In a series of earlier papers we studied the tidal response of an early type star in a close binary system (e.g. a massive X-ray 
binary), ignoring stellar rotation, but taking into account the effects of internal nuclear evolution (Savonije & Papaloizou 
1983, 1984: SP83,SP84 and Papaloizou & Savonije 1985: PS85). Recently, we extended these studies by taking into account 
the effects of the Coriolis force in uniformly rotating stars (Savonije, Papaloizou & Alberts 1995: SPA95). We found that, 
as expected, in a rotating star resonances between the I = 2 tide and internal modes are possible not only for I = 2, but 
also with modes with predominantly Z = 4 or Z = 6 etc. due to the mode coupling by the Coriolis force. Once the relative 
forcing frequency a was smaller than twice the stellar rotation speed, however, the numerical results were swamped by short 
wavelength inertial waves. We could only study the (positive frequency range of the) inertial regime by adopting the so called 
‘traditional approximation’ in which the 0-part of the stellar angular speed Qa is ignored. It appeared that rotationally modified 
g-modes continued to be excited in this frequency range but that for a < fig, no reliable results could be obtained. Therefore 
a new implicit 2D scheme was developed that is numerically more stable and that can treat three-level difference equations, 
enabling viscosity to be introduced in the equations of motion. The viscosity can then be used to damp the shortest wavelength 
inertial waves. We have used this new scheme to explore the inertial regime, introducing a much finer grid resolution in the 
radial direction. For free oscillations, the low frequency inertial regime was studied by e.g. Berthomieu et al. (1978), Saio 
(1982) and by Lee & Saio (1987, 1989). However, these authors did not take the Coriolis force fully into account, but used 
the ‘traditional approximation’, or truncated (to the first few terms) the infinite expansion in spherical harmonics required 
to describe the oscillations in rotating stars. We take the Coriolis force fully into account, by not using spherical harmonics 
expansions, but solving the full 2D problem in r and 9. Also we consider (tidally) forced oscillations. The next section (2) 
summarizes our basic equations, followed by a brief section (3) on the stellar input model. Then, in section 4, we present our 
numerical results. In the appendix we list the finite difference approximation to the partial differential equations describing 
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non-radial oscillations that give stable numerical results and a short description of the new implicit 2D elimination scheme 
used to solve these equations. 


2 BASIC EQUATIONS 

We consider a uniformly rotating, early type secondary star with mass Mg and radius Rs in a close binary with circular orbit 
with angular velocity lo and orbital separation D. We assume the secondary’s angular velocity of rotation Dg to be much 
smaller than its break-up speed, i.e. (flg/flc)^ 1 , with = GM^/R^, so that effects of centrifugal distortion (oc flg^) 
may be neglected in first approximation. We wish to study the response of this uniformly rotating star to a perturbing 
time-dependent tidal force. The Coriolis acceleration is proportional to fig and we consider its effect on the tidally induced 
motions in the star. We use spherical coordinates with origin at the secondary’s centre, whereby 6 — 0 corresponds 

to its rotation axis which we assume to be parallel to the orbital angular momentum vector. We take the coordinates to be 
non-rotating. 

As is well known, in a non-rotating star the solutions of the linearized non-radial stellar oscillation equations can be 
expressed in terms of spherical harmonics, i.e. the spatial part of each mode can be fully separated into r-, 9- and (^s-factors 
(e.g. Ledoux & Walraven 1958) 

U{r,9,ip) = u{r)PP{cos 9)e’^'^, 

where represents the associated Legendre polynomials for I and |m|. 

The introduction of the Coriolis force, however, destroys the full separability of the oscillation equations, it only being 
retained for the </ 5 -dependence. It turns out (e.g. Berthomieu, Gonczi, Graff, Provost & Rocca 1978) that two independent sets 
of approximately spheroidal oscillation modes exist: modes in which the density perturbation is even with respect to reflection 
in the equatorial plane, which have I — \m\ even valued, and modes with odd symmetry for the density, having I — \m\ odd 
valued. In addition, for each I, there is a set of approximately toroidal modes (e.g. Papaloizou & Pringle 1978) which couple 
with the spheroidal modes of Z ± 1 . 

Let us denote perturbed Eulerian quantities like pressure P', density p', temperature T' and energy flux F' with a prime. 
The linearized hydrodynamic equations governing the non-adiabatic response of the uniformly rotating star to the perturbing 
potential may then be written 



where are the unit vectors of our spherical coordinate system, k is the unit vector along the rotation axis, v' denotes the 
velocity perturbation, k the opacity of stellar material and S its specific entropy. These perturbation equations represent, 
respectively, conservation of momentum, conservation of mass and conservation of energy, while the last equation describes 
the radiative diffusion of the perturbed energy flux. For simplicity we adopt the Cowling (1941) approximation, i.e. we neglect 
perturbations to the gravitational potential caused by the secondary’s distortion. We also neglect perturbations of the nuclear 
energy sources and of convection. 

The dominant tidal term of the primary’s perturbing potential is given by the real part of 

$T(r, 0, p, t) = fp Pi(p) (5) 

where p = cos 9, PUp) is the associated Legendre polynomial for I = |m| = 2 and / = with Mp the companion’s mass. 

Adopting the same azimuthal m = 2 symmetry and time dependence for the perturbed quantities as the forcing potential, 
the perturbed velocity v' = 2i(a; — Og)^, where ^ is the displacement vector. The perturbations can be written as e.g. 

Gir,9,ip,t) =^r{r, 9) 

where is the radial component of the displacement vector and is assumed complex to describe the phase shift A 93 

with respect to the forcing potential (^). 

For numerical reasons it proved better not to eliminate from the equations, as was done in SPA95, and retain the (p- 
equation of motion. However, we do eliminate with help of the (^-component of the radiative transport equation. Writing 
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for simplicity from now on for ^r{r,0) etc., while dividing out the factor , equations (0)-(^) transform into the 

equations given below. 

The radial equation of motion becomes 


pcT^r + [2i(jf2s sin 6p] 






6fpr (sin 6)^ = 0 


( 6 ) 


where a = 2uj = 2(a; — rig) is the forcing frequency felt by a mass element in the uniformly rotating secondary, while 
Xp = gif InP and Xt = gifr follow from the equation of state. 

The 0-equation of motion becomes 


pa% + l2iansCOsdp]^^- (^) ^ ^ (y) - 6/pr sin0 cos 0 = 0 

The (p-equation of motion becomes 


per 


picrfla sin 0p] 


[2i(Tns cos 0p] + 


"imPxp ' 

rsin0 


^ f imPxT ' 
\P ) L sin 0 


[t 


+ Simfpr sin 0 = 0 


( 7 ) 

( 8 ) 


The perturbed equation of continuity is transformed into 

p' 


P . 1 1 

- 1 - — VT Pir ) + ■ 

p r^p or ^ ' 


^ [sineio)-= ^ 


rsin0 90rsin0'’'^ 
By applying the thermodynamic relation 


( 9 ) 


5S = S' +i.VS 


P 1 (^P Y ^ 

pT Ts - 1 V P % 


where the symbol 5 denotes a Lagrangian perturbation and Tj the adiabatic exponents of Chandrasekhar, the perturbed 
energy equation can be transformed into 


rdlnP 
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( 10 ) 


is a local characteristic radiative diffusion length in the star, with F the unperturbed (radial) energy flux. When p —> 0 
the response becomes adiabatic and shows (outside resonance) no phase shift with respect to the companion. However, even 
for ‘high’ frequencies a/Qc ~ 1 the diffusion length p becomes comparable to the scale height when the stellar surface is 
approached, so that the imaginary part of the surface response becomes comparable to the real part. 

The perturbed radial energy flux is given by 







where Kp = gjf^ InK and kt = g-§p^ Infr. Finally the 0-component of the perturbed flux follows as 


( 11 ) 



where p = cos 0. 


2.1 Boundary conditions 


The differential equations are supplemented by the following boundary conditions: at the stellar centre we require and Ft 
to vanish, while at the stellar surface, we require the Lagrangian pressure perturbations to vanish 


P' dlnP^ „ 

P dr “ 

and the temperature and flux perturbations to fulfill Stefan-Boltzmann’s law 


(13) 


F' 

Y 



dlnT\ 

dr 7 


Cr + 4 


{ n. 

[t 


( 14 ) 


© 0000 RAS, MNRAS 000, 000-000 

































4 G.J. Savonije and J.C. Papaloizou 

Furthermore, and Fg must vanish on the rotation axis while in view of the symmetry of the tidal force we adopt mirror 
symmetry about the equatorial plane, i.e., for 9 = tt/2 we also require ^g and Fg to vanish. 


2.2 Viscous damping of short wavelength response in convective regions 


In order to damp the short wavelength (of order gridsize) response to the tidal forcing in the inertial regime in the non- 
stratified (convective) regions of the star, viscous terms are added to the equations of motion. We use a simple mixing length 
prescription for both energy transport and turbulent viscosity in convective regions. We adopt a (constant) reduction factor e 
for the turbulent viscosity to roughly take into account the mismatch between the forcing time and the convective timescale 
in the core. The turbulent viscosity can thus be expressed as 


C = e Ffp 



(15) 


where F is the local energy flux and Hp the local pressure scale height. For the turbulent viscosity in the convective core a 
reduction factor e = 10“^ is adopted (adequate to damp the grid oscillations for all applied forcing frequencies), while in the 
convective outer shell a slightly larger value e = 3 x 10“^ had to be taken. At the boundary of the convective regions we let the 


viscosity C, drop to negligible values as exp 




, where d is the distance from the boundary. Fig. hi shows the adopted 


(( distribution through the star. A significant steeper decline to zero viscosity at the convective core boundary would not 
prevent the short wavelength inertial response to penetrate into the weakly stratified regions adjacent to the core. Adopting 
a much steeper decline (decline on scale O.lHp) changed the tidal torque by less than a few percent. In the convective surface 
shell near i ~ 1600, the ratio of the forcing time to the local viscous time (where Ar is the local gridsize) becomes 

larger than unity. But since the mass in these layers is negligible this has no effect on the tidal torque, which is predominantly 
applied in the radiative region adjacent to the convective core. 

With this prescription for the viscosity in the convective regions the grid oscillations are adequately damped while the 
global numerical results, like the tidal torque, appear insensitive to the adopted viscosity. Only if ^ is taken some orders of 
magnitude larger do the global numerical results like the net tidal torque change significantly. 

The (simplified) viscous terms that are added to the equations of motion are expressed in such a way that the PDF’s 
(^)-(^) retain their separability in r and 9 when fig ^ 0. 

The radial equation of motion (^ is extended with the following viscous terms 





Ofl 




(16) 


while to the 0-equation of motion the following terms are added 





5- 





(17) 


3 THE UNPERTURBED STELLAR MODEL 

A recent version (Pols, Tout, Eggleton & Han 1995) of the stellar evolution code developed by Eggleton (1972) was used to 
construct the unperturbed stellar input model for the secondary. The model represents a spherical zero-age main-sequence star 
of 20 Mq with chemical composition X — 0.70 and Z — 0.02. The model was constructed with the OPAL opacities (Iglesias 
& Rogers 1993). With these opacities the 20 M© ZAMS star appears to have shallow convection zones in the envelope which 
are absent when the Cox-Stewart opacities are used. The stellar radius equals Rb = 4.993 x 10^^ cm, while the stellar moment 
of inertia 7s = 7.21 x 10®® g cm^. The break-up angular speed equals He = 1-46 x 10““* s“^. In the convective core, which has 
a mass of about 10 Mq, the imaginary Brunt—Vaisalla frequency A'^ given by 

^ nsi 

p dr ^ p dr 'yP dr J 

1 _q 1 

reaches a maximum value of about |A| = = 5 x 10 ^ Hs, where \A\'^ is given in units of He, while in the outer convective 

shell 1A| 2 ~ 5, see Fig. ^ This corresponds to a super-adiabaticity V “ Vad — 10~® and 0.07 in the core and envelope shell, 
respectively. 

It can be seen that the radial grid has about 500 zones in the convective core, about 600 zones in the adjacent radiative 
‘plateau region’ where the eigenfunctions of g- and r-modes have most of their radial nodes, and also several hundred zones 
in the rarified outer layers with convective shells. 
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stellar input model 



Figure 1. Characteristics of the 20 Mq stellar input model: the (adapted) Brunt—Vaisalla frequency i/bv = signyty^jyt| in units of the 
stellar break-up speed flc (vertical scale on the right), radius r in cm and turbulent viscosity in cgs units (vertical scale on the left) as 
a function of radial mesh number. 


3.1 Non-adiabatic surface layer (NASL) 

A crude random walk estimate for the time required for heat to diffuse to the stellar surface from a given point r in the stellar 
envelope is given by 

Pk{Rs - r)'^f3 
“ c(l - P) 

where k is the opacity, c the velocity of light, and /3 the ratio of gas to total pressure. By comparing this timescale with 
the forcing time r/ = 27r/a, we find that non-adiabatic effects become important (rth — Tf) roughly from radial zone 1300 
(for Z<7 = O.Olflc) or 1400 (for oJ = 0.20flc) outwards. The convective shells are thus located in this thin non-adiabatic layer, 
see Fig. so that the (adiabatic) Brunt—Vaisalla frequency | A| 2 greatly overestimates the convective growth rate there 
{rpf ~ 3flc in the innermost convective shell). Also the thermal stratification in the radiative part of the non-adiabatic layer 
is much weaker than suggested by the large values |A|^ ~ 5flc or even higher, attained in that region. 


4 NUMERICAL RESULTS 
4.1 Introduction 


We extend the calculations of SPA95 into the low-frequency inertial regime, i.e. for forcing frequencies |ct| < 2fls. In SPA95 the 
numerical results in this frequency regime were swamped by short wavelength inertial waves in the convective core, which also 
propagated into the weakly stratified radiative regions adjacent to the core, unless we adopted the ‘traditional approximation’ 
in which the 9 component of the stellar angular velocity fig is neglected. And even then, we could not obtain reliable results 
for oJ < 0.1, where oJ = is the tidal pattern speed in the corotating frame expressed in units of flc. From now on we shall 

express all frequencies in units of flc = {GMa/Rs^)^. _ 

Equations (^- (^) have been approximated by the set of finite difference equations (FDE’s) (M given in the appendix. 
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The FDE’s with appropriate boundary conditions are then solved by a newly developed 2D implicit elimination scheme briefly 
explained in appendix B. With this new three-level finite difference scheme we can effectively damp the very short (gridsize) 
inertial waves by introducing a modest (turbulent) viscosity in the equations of motion, see equations 


4.2 Fourier-Legendre coefficients 

For each forcing frequency o’ the calculated stellar response can be expanded in Fourier-Legendre (FL)series with complex 
coefficients Ci^{r), e.g. 


G{r,n) = ^cr{r)Pr{n) 

1=0 


with 




'(u)du. 


By plotting the real and imaginary parts of (r) for the various perturbed variables versus r for m = 2 and different Lvalues, 
we can compare the strength of the different /-components in the response. 


4.3 Tidal response in non-rotating stars 


To compare our results for (uniformly) rotating stars with those of non-rotating stars we have recalculated (for w < 0.1) the 
tidal response for a non-rotating star with the same 2-D scheme and with the same viscosity. For forcing frequencies in the 
range considered here the amplitude of the displacement vector ^ globally follows that of the equilibrium tide, except that in 
the subsurface layers {i > 1500), a peak occurs which grows with decreasing tidal frequency o = 2a;. For u) = 0.025 the radial 
displacement peaks at more than twice the value of the equilibrium tide near i = 1500, while decreases to the equilibrium 
value near the stellar surface. The 6'-displacement is generally an order of magnitude larger than near the surface. The mass 
in the NASL is only small, so that the effect on the total torque is generally insignificant. The likely rea son for this strong 
sub-surface response is that for low forcing frequencies the radiative losses in the surface layers (section 3.1) get so strong 
that the characteristic timescale for convective oscillatory motions becomes comparable to the forcing period. Indeed, when a 
zero value for the radiative diffusion length r) (adiabatic calculation) is adopted the strong sub-surface response disappears. 


4.4 Tidal response in sub-synchronous stars (fig < a;) 


For sub-synchronously rotating (ZU > 0) early type stars rotationally modified g-modes are excited in the radiative envelope 
for all frequencies in the inertial regime, down to near co-rotation. The g-modes have a predominant I = 2 component, with 
I = 4 and 6 components having in general a factor 4-10 weaker amplitudes in this frequency range. If the relative forcing 
frequency a = 2Zj = 2(a; — Dg) in the frame corotating with the secondary falls below 2IIg, inertial waves are excited in the 
weakly or non-stratified layers. Once ZF < fig the first signs of inertial waves are observed, not only in the convective core, but 


even stronger in the non-adiabatic surface layer (NASL) containing the convective shells (section 3.1). Without viscosity the 


tidal response in the convective regions gets overwhelmed by short lengthscale (gridsize) variations, as discussed in SPA95. 


The excitation of inertial modes in the core is discussed in section 1.7 


As discussed above for non-rotating stars, the strong radiative damping in the NASL gives rise to large amplitude 
oscillatory convective motion in the sub-surface layers. With rotation the non-adiabatic response in the sub-surface layers 
shows even larger ^-displacements, while the response is more complex with high /-components due to inertial effects. 

Fig. I shows an example of the ^-dependence of the tidal response in the sub-surface layers of a star rotating with Dg = 0.1 
for a forcing frequency ZU = 0.076. At i = 1400 (not shown) the 9 dependence is still typically that of the radiative region 
outside the convective core and consists predominantly of a simple 1 — 2 (^r) or / = 3 (^e) component. The displacement ^e/Ra 
attains a maximum of about 8 x 10~'* at i = 1400. However, slightly further out (from about i = 1500) inertial waves are 
excited which show up as short wavelength oscillations. The amplitude of the displacement vector ^ shows a peak at around 
cos6' ~ ZZJ/Dg with an oscillatory decay towards the rotation axis. Note that g-modes become evanescent when the oscillation 
frequency is smaller than twice the average radial component of the stellar angular velocity 2ng cos 6 (e.g. Eckart 1960), so 
that the peak amplitude of the tidal response appears to coincide roughly with the turning point of those modes. Inertial 
waves can propagate from pole to equator. 
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0/n 0/n 


Figure 2.u3 = 0.076 and Hs = 0.1: the 0-dependence of the radial- and 0-component of the displacement vector for radial zones in the 
range i = 1650 — 1750, corresponding to r/Rs = 0.992 — 0.998, i.e. in between the two convective shells near the surface. The left hand 
panels give the real parts, the right hand panels the imaginary parts. The shortest lengthscales were removed by viscous damping. 


4.5 Tidal response in super-synchronous stars (fig > w) 

When we lower the forcing frequency - through corotation to negative values - the tidal response becomes qualitatively different. 
Now the excited tidal wave becomes retrograde as the star is spinning faster than the orbital revolution of its companion. For 
small |IU| there are again many radial nodes in the radiative envelope, as there are for small positive frequencies, but now the 
displacements show a strong toroidal component. When |ZU| is increased a rapid succession of strong resonances occurs, these 
are resonances with free r-modes, e.g. Papaloizou and Pringle 1978, Provost, Berthomieu & Rocca (1981), Rocca (1982) and 
SPA95, analogous to Rossby-waves in the earth’s atmosphere (e.g. Pedlosky 1977). Because the excited r-modes occupy only 
a small frequency range approximately given by — < uJ < 0 with I = 3 and m = 2, strong resonances with low order 

(with few radial nodes in the envelope) r-modes are possible for relatively small |Z<7|, i.e. close to corotation. 

Hence, in contrast to early type sw6-synchronous stars approaching corotation, which tend to spin up towards corotation 
as a result of weak resonances with high-order, strongly damped g-modes, slightly super-synchronous stars tend to spin down 
towards corotation much more efficiently due to strong resonances with weakly damped r-modes (see Figs. 

Figs. ^ shows the radial distribution of the dominant FL-coefficients of the radial and 0-displacement for the strong 
resonance at ZiJ = —0.0166 with the lowest order r-mode that could be found for a star rotating at fls = 0.1. In the radiative 
envelope the radial displacement is typically three orders of magnitude smaller than the 0-displacement which shows the 
(quasi- (toroidal character of this mode. 

For this low order r-mode resonance, |^| attains values large compared to the equilibrium tide (Fig. ^ not only in the 
radiative envelope, but also in the convective core. The response in the core is oscillatory, both in r and 0 (Fig. ^ indicating 
the excitation of an inertial core mode. Fig. ^ shows the 0-dependence of the radial and 0 displacement in the convective shell 
region. 

When Zu is made more negative than ZiJ ~ the imaginary parts of the perturbed quantities decrease substantially. 

The imaginary part of p' shows weak oscillations in the core and envelope, with an evanescent zone in between (from about 
i=700 to i=1100). At the boundary of the convective core a (relatively) strong spike occurs. When tU is made more negative, a 
very short wavelength response appears just outside this spike, which starts to extend outwards, grow in amplitude, and fill up 
the ‘evanescent’ region. At Zj = —0.0375 (for fig = 0.1) some 80 radial nodes occur in the radiative envelope. It appears that. 
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RADIUS RADIUS 

Figure 3. ZZJ = —0.0166 and fig = 0.1: corresponding to a stellar rotation period of Ps = 5.0 days and orbital period ^ 6 days. Plotted 
are the dominant Fourier-Legendre coefficients versus radius for the radial and ^-displacement for the lowest order resonance with an 
r-mode. The various curves correspond to the consecutive ^-values shown inside each panel (from left ^-value to right ^-value: full line, 
dashed line, dot-dash, dotted and dash-dot-dot-dot). The real number C(l) given after each value is the maximum value attained over 
the star for that particular FL-coefficient. In the panel the dash-dot-dot-dot line gives the corresponding amplitude of the equilibrium 
tide. The surface layers {i > 1300) have been cut away in order to see the smaller amplitude structure of Im(^r) in the interior. 


by decreasing the forcing frequency, a transition has been made from the r-modes to the rotationally modified, retrograde 
g-modes. For more negative frequencies the number of radial nodes decreases and we proceed through resonances with lower 
and lower order g-modes, until the limit ZJ = —fig is approached when u; ^ 0. 


© 0000 RAS, MNRAS 000, 000-000 














ivective 


depend 
by vise 














10 G.J. Savonije and J.C. Papaloizou 


n^= 0.10 u-Q 3 = 0.0185 to= 0.411E+02yr Torque= 0.127E+37cgs 


ro 

O 

X 

CM 


O 

a> 

(D 


O 


o 

X 

CM 

I 


0.2 0.4 0.6 0.8 

r/R 


1 



Figure 6. The torque integral ( p^ ) integrated from 0 to r as a function of r for ZJ = 0.0185 and Qb = 0.1. This 
corresponds to a strong resonance with an inertial mode in the core, causing a broad dip in the to curve, see Fig. 
plotted torque integral are given in units of GM^/Rg. 


forcing frequency 
The values of the 


4.6 Tidal torque and stellar spin-up/down time 


The radiative damping in the surface layers of tidally excited non-radial oscillation modes introduces a phase shift between 
the mode and the forcing potential. This phase shift gives rise to a net tidal torque on the secondary 


Te* 



Im{p') Ppi(n) drd/r 


(19) 


When the secondary rotates slower than synchronously (fig < w) the net torque is positive and give rise to a spin up, 
while for super-synchronous stars (fig > ui) the tidal wave is retrograde in the frame co-rotating with the secondary and yield 
a negative tidal torque, leading to spin down. 

The secondary’s spin-up/down time-scale tap follows as 

, — 1 "Pext 

tap ~ r — 

where la is the secondary’s moment of inertia. Let us now define, as in our earlier studies, the ‘intrinsic’ spin-up time to of 
the secondary as 

to(aJ)= (^)'(§)%sp(cF). (20) 

in which we have divided out the dependency of the I = m = 2 forcing on the specific binary configuration. 

By varying the forcing frequency o = 2lJ and calculating the tidal torque according to equation (^^ we can determine the 
‘intrinsic’ tidal spin-up time to (equation §) as a function of to. 

Fig. 0 shows the calculated to distribution for a non-rotating star, whereby the numbers along the curve indicate the 
number of radial nodes of the tidal oscillations. The dotted curve corresponds to a low-frequency approximation (PS85): 


to{uj) 


flcVMgJ \Ra) Vfic/ 


( 21 ) 


where k is the star’s gyration radius {la = kMsRa^) and Rcc, Mcc are the radius and mass of the convective core. This 
expression is based on a simple model which assumes a purely adiabatic star with a thin very non-adiabatic outer shell, so 
that all outward propagating gravity waves are absorbed by radiative damping in the outer shell. Expression (^) is similar to 
that obtained by Zahn (1977) and is valid only for non-rotating unevolved massive stars in the limit of low forcing frequencies. 
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Figure 7. The scaled spin-up or spin-down time to versus lo for a non-rotating 20 Mq star. For comparison the low frequency limit 
is also plotted (dotted curve). The number of radial nodes is indicated along the curve. 


21D 


In Papaloizou & Savonije 1997 (PS97) the asymptotic low frequency expression for to is generalized to include the effects of 
stellar rotation. 

Figs, ^and § show the calculated intrinsic spin-up or down time to(tn) versus tidal frequency for a 20 Mq star rotating at 
a rate Qa = 0.1 and IIs = 0.2, respectively. For low frequency g-modes, which have many radial nodes, a strong cancellation 
effect occurs whereby the torque integral, when integrated over the envelope, changes sign many times. This can be seen in 
Fig. § where the torque integral dii), integrated from r = 0 to r, is plotted versus r for the response at w = 0.0185. At this 
frequency a strong resonance occurs with an inertial mode in the core, see next section. Meshpoints 1500-1800 (the convective 
shell region) correspond to r/R > 0.98. Note that for small positive values of the forcing frequency the wavelength of g-modes 
in the envelope becomes too small (> 50 radial nodes) for adequate resolution on the grid. 

For rotating stars which can have negative forcing frequencies when the star spins faster than the orbital revolution of 
its companion, the to curve has an extra branch for ZU < 0. Rotational effects make the to curve qualitatively different, not 
only for the new branch with ZZ7 < 0, but also for ZZ7 > 0. The Coriolis force couples modes with even I— values, so that I = 2 
forcing indirectly excites Z = 4,6 .. as well. It also gives rise to a frequency shift of the non-radial oscillations. Comparing 
Figs. 0 and 1 ^ shows that, as expected, a shift to lower eigen-frequencies has occurred in the rotating star (e.g. Ledoux & 
Walraven 1958, SPA95). The response for negative forcing frequencies is distinctly different from that for positive frequencies. 
For small negative forcing frequencies resonances with free rotationally controlled (quasi-)toroidal r-modes show up as strong 
narrow dips in the to curve. These strong resonances can cause a rapid spin-down of slightly super-synchronous stars. For 
more negative forcing frequencies once again resonances with (retrograde, rotationally modified) g-modes become possible. 
The more negative forcing frequencies correspond to lower order resonances (with fewer radial nodes) that are damped only 
weakly, with consequently smaller to values. 
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Figure 8. The scaled spin-up or spin-down time to versus uJ for a uniformly rotating 20 Mq star at a rate Qs = 0.1. The number of 
radial nodes in the envelope is indicated along the curve. For comparison the low frequency limit (pl|), whereby we have simply replaced 
uj by oJ, is also plotted (dotted curve). 





Figure 9. The scaled spin-up or spin-down time to versus uJ for a uniformly rotating 20 Mq star at a rate fls = 0.2. The number of 
radial nodes in the envelope is indicated along the curve. For comparison the low frequency limit (^), whereby we have simply replaced 
Lj by lJ, is also plotted (dotted curve). 
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Figure 10. The dominant FL-coefBcients for the radial and 9 displacement versus radius for the tidal response at ZU = 0.0185 and 
f!a = 0.1 (similar to figure |^. For this forcing frequency the tide is resonant with an inertial mode in the convective core. The convective 
core boundary is at i ~ 460. The non-adiabatic surface layer, where the amplitude becomes large, is not shown. 


4.7 Excitation of inertial modes in the convective core 

Although the inertial modes have a dense spectrum their amplitude under tidal forcing remains limited. This is probably 
partly due to coupling with envelope modes which keeps the amplitude down because of leaking. However, the long scale tidal 
force is not efficient in exciting inertial modes. For a constant density and pressure core it can be shown (PS97) that only 
two potential inertial mode resonances are excited by the the tide, for tj/fls = (±Vl05 — 7)/28. For fig = 0.1 these roots 
correspond to ZF = 0.012 and ZF = —0.062. For a core with a density and pressure gradient, like the one studied in this paper, 
inertial modes can be excited throughout the inertial range, but generally with an amplitude smaller than the response for 
the two above mentioned resonance frequencies. 

Note that in Figs. ^ and which show the to distributions obtained by the detailed numerical calculations for the more 
realistic stellar model considered here, two broad resonance dips are visible corresponding roughly with the above mentioned 
two resonant frequencies. In Fig. ^ we have plotted, for the positive resonance frequency ZF = 0.0185 (fig = 0.1), the FL- 
coefficients of the radial and 9 displacement versus radial mesh number. It can be seen that the imaginary parts of the resonant 
inertial response in the core attain amplitudes large compared to those of the (high order) g-mode in the radiative envelope. 
Apparently the inertial mode in the core excites g-modes in the envelope to a level higher than normal, yielding the broad 
dip in the to curve. 

At ZF ~ 0.007 another strong resonance occurs with an inertial mode in the core. However, the resolution in the envelope 
is too poor to determine to for such low frequency response. 


© 0000 RAS, MNRAS 000, 000-000 







14 


G.J. Savonije and J.C. Papaloizou 



Figure 11. The tidal spin-up or spin-down timescale tsp as a function of Rb/D, where D is the orbital separation, for a binary system 
consisting of two 20 Mq unevolved stars. The full fat curve corresponds to a stellar spin rate fig = 0.2, the dotted curve to Qb = 0.1 and 
the thin full line to a non-rotating star. The spin up-time curve for fig = 0.2 correspond to larger values of Rb/D than shown here. 


4.8 Comparison with observations 

Giuricin, Mardirossian and Mezetti (1984 a,b,c) inspected a large sample of early type binaries and studied their orbital 
eccentricity and synchronism between orbital revolution and stellar spin rate. Their conclusion is that Zahn’s asymptotic 
expression (1977) for tidal effects can explain the majority of observed eccentricities. There remains, however, a not insignificant 
fraction of relatively wide early type binaries with nearly circular orbits that still require an explanation. Furthermore, these 
authors note a pronounced tendency towards synchronisation up to relative radii Rb/D ~ 0.05, while Zahn’s asymptotic 
theory predicts a limiting separation for this corresponding to Rs/D ~ 0.15. 

To synchronise a wide (circular) system with equal masses and Rb/D = 0.05 during the main sequence stage of duration 
tms of its early type components requires a time averaged ‘intrinsic’ spin-down time < to >— (0.05)® tms, which for 
tms — 8 X 10® year corresponds to < to >— 0.1 year. Rb/D = 0.10 would require < to >— 10 years. That demands a strong 
tidal response, which is difficult to attain with radiative damping of dynamical tides, unless the tide is in continuous resonant 
interaction with a weakly damped mode. But note that in a wide binary with Rb/D = 0.05 and mass ratio unity (which 
corresponds to Pi, ~ 30 days) an early type star rotating with a period of 2.5 days (fig = 0.2) encounters strong resonances 
with weakly damped, rotationally modified, retrograde g-modes for which to — 1 year (Fig. ^). Hence it is expected to undergo 
significant spin down. 

To visualize our numerical results we have plotted in Fig. the spin-up or spin-down time as a function of Rb/D for a 
circular binary system consisting of two unevolved 2 OM 0 stars, by scaling the general results to{uj) to this particular binary 
configuration, adopting fls = 0.2 (full fat curve), fis =0.1 (dotted curve) and fig = 0 (thin full line). The spin-up times for 
fig = 0.2 are not shown, they correspond to Rb/D values larger than given in the figure. 

It can be seen that with rotation the spin down time in wide binaries tend to be shorter, i.e. tidal effects can be significant 
even in binaries with R/D ~ 0.06. ft should be taken into account that the star is subject to internal (nuclear) evolution 
causing the core to shrink and the envelope to expand during the core hydrogen burning stage. The envelope expansion causes 
the star to spin down by about a factor 2 during the main sequence stage (this factor would be larger when the envelope 
decouples from the more rapidly spinning core). The core contraction, on the other hand, gives rise to increased gravity and 
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a shift of the g-mode frequencies to higher (absolute) values. Thus the internal evolution has signihcant effects on the tidal 
evolution of the binary (SP84). When a mode frequency shifts to higher frequency the radiative damping decreases and the 
corresponding to goes down, making it easier to synchronise the star. 

What could further ease synchronisation is that the tides might induce differential rotation whereby the surface layers 
get synchronised and decoupled (to some extent) from the more rapidly spinning interior. Then apparent synchronisation of 
a main sequence star in quite wide binaries is relatively easy. Goldreich & Nicholson (1989) suggested that this is actually to 
be expected because the tidal despinning takes place where the g-modes are dissipated, i.e. initially in the surface layers of 
the star. When the surface layers have been synchronised, the gravity waves (excited near the core boundary) can no longer 
penetrate to the surface. When gravity waves approach the corotating surface layers their wavelength keeps diminishing until 
they are fully absorbed. The absorbed negative angular momentum carried by the gravity waves then leads to further spin 
down of the layers immediately beneath the corotating layers. Adopting this picture, the more rapidly spinning interior keeps 
its initial rotation rate with a slowly rotating, synchronized shell on top of it. Shear instabilities or other processes of angular 
momentum redistribution are expected to smooth this discontinuous angular velocity prohle to some extent. Assuming that 
the to curve for the original uniform rotation rate remains approximately applicable, we deduce from Fig. |^that for fig = 0.2 
the spin down time tsp (required to reach corotation) is on average roughly 10 times longer than the star’s main sequence 
lifetime for systems with 0.15 > R/D > 0.06. Hereby we take into account that structural changes of the star due to stellar 
evolution shifts its modal spectrum, so that the star cannot remain stuck on a bump in the to curve. This would correspond 
to synchronisation of the outer several percent of the stellar mass, i.e. synchronisation down to a depth oi r/Rs — 0.6 during 
the main sequence phase. It suggests that tides in rotating stars may cause ‘apparent’ synchronisation in quite wide binary 
systems. 


5 DISCUSSION 

We have studied the fully non-adiabatic tidal response of a slowly and uniformly rotating 20 Mq ZAMS star to an 1 = m = 2 
perturbing potential in the low frequency, inertial regime. We solved the partial differential equations governing the response 
using a 2D-implicit hnite difference scheme. In the work described here, the Coriolis force was included but the centrifugal 
force was neglected. This enabled us to work with spherically symmetric equilibrium models. 

Our results show that, as expected, rotational effects can signihcantly influence the tidal response. We have extended the 
study of tidal effects to super-synchronous stars in which retrograde r-modes and rotationally modihed retrograde g-modes 
are excited. The r-modes give rise to efficient spin-down in slightly super-synchronous stars, while the retrograde g-modes can 
give efficient despinning of highly super-synchronous early type stars in wide binaries. For all frequencies in the inertial regime 
inertial waves are excited in the core and convective surface shells. Consistent with analytic results derived for a simplihed 
model with a constant density core (Papaloizou and Savonije 1997) the detailed numerical results show the existence of 
frequency intervals for which the tidal response is significantly enhanced due to resonance with an inertial mode in the core. 
The calculated broad resonance features indicate that the resonant inertial mode in the (almost adiabatic) core is damped 
by exciting gravity modes in the envelope to a level higher than normal. Assuming that the tides induce differential rotation 
by synchronising the outer parts of the early type star, while leaving its interior roughly unaltered, synchronisation may be 
achieved in a substantial part of the envelope during the main sequence phase, even in quite wide binary systems. This could 
explain the observed approximate synchronous rotation in some wide binary systems with R/D ~ 0.05. 
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Figure Al. Schematic picture of the elemental gridcel for the three-level finite difference scheme used. The lines j =constant all intersect 
at the stellar centre. The symbols ©, □ and O denote the location where the different variables are defined and the locus where the 7 
difference equations are centred. 

© corresponds to and centre of equations (A1,A6) 

□ corresponds to and Fq ; centre of equations (A2,A7) 

O corresponds to p' and T' ; centre of equations (A3,A4 and A5) 


APPENDIX A: THE FINITE DIFFERENCE EQUATIONS (FDE’S) 

We solve, as function of r and 6^ the set of seven linearized partial differential equations, together with boundary 

conditions, for the seven complex quantities ^r, ^ 9 , ^ip, rewriting them as a set of finite 

difference equations on a polar (r, 6) grid in a meridional plane of the star. Because the forcing has symmetry about the 
equatorial plane {6 = 7r/2) we need only consider the hrst quadrant: Q < 6 < -k 12. The index i = l,Nr runs along radius (from 
r = 0 to r = Rs) and the index j = 1, Ng runs in the 0-direction (from 0 = 0 to 0 = 7r/2). 

Note that for stable numerical results care is needed to construct the FDE’s. The set of FDE’s given below yielded 
satisfactory numerical results. These FDE’s approximate the PDF’s in such a way that for fls ^ 0 they too become separable 
in r and 0. For numerical stability we use a staggered mesh, so that the various perturbed quantities to be solved for are 
not all defined at the same location. See Fig. Al for a schematic picture of an elemental gridcel for the three-level difference 
scheme used in the calculations. 

The adopted radial grid is not equidistant. The 0-grid is for 0 > ^ equidistant in /i = cos 0 with A/r = while it 

is for 9 < J equidistant {Av = i^j — Vj-i) in v = sin0. With this 0-grid we have for Ng = 64 adequate resolution, especially 
both near the rotation axis and near the stellar equator. At the transition point j = Nh (corresponding to 9 = ^) v = 
However, because of the staggered mesh, for 2"“* order derivatives expressions for a non-equidistant (0-)grid should be taken 
at that point. 

The radial equation of motion (& is approximated by 
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The 0-equation of motion (^) is approximated by 
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Finally, the energy diffusion in the direction ( |l^ ) is approximated by 
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where A was defined in equation (A2). 
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A1 The viscous terms 

The viscous terms ( 00 ) added to the equations of motion contain second order derivatives and must therefore be defined 
on three levels of the grid. Again the finite differences are constructed in such a way that the expressions retain separability 
in r and 6. To that end the derivatives with respect to in equations (00 are for j < Nh expressed as 

d^r 


_ _ _ tLT~ _ 
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where for ffg = 0 we have ie = ^ (x u. With at = iaRa{i/r‘^)i, the radial derivatives to be added to the equation of radial 
motion & become 
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where Ai = n+i — ri, As = n — ri_i and A 3 = n+i — n-i. 

The 0-derivatives to be added to the equation of radial motion become for j < Nh 
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where Ai = v-: 1 — , As = 1 — Vi_3 and A 3 = v-, 1 — u- 3 . Note that at the transition point j = Nh the grid is 
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not equidistant (neither in jj, nor i/) and that it is essential to use the expressions for non-equidistant mesh given here. For 


j > Nh these derivatives become 
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where now Ai = Uj+l ~ Pj-ii ^2 = Mj-l “ '^3 ~ ^j+- ~ similarly for the 0-equation of motion. 


APPENDIX B: NUMERICAL SOLUTION TECHNIQUE 


Let us introduce the following definition for the perturbation vector Xij' at each gridpoint centre consistent with Fig. AI 
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where F, p and T respectively denote the local unperturbed energy-flux, density and temperature and Rs the stellar radius. 

Then the finite difference equations (A1-A7) corresponding to the 7 linearized partial differential equations (^) —(|l2|) 
extended with the viscous terms can for each step {i,j) be expressed as a single matrix equation 


AiXi_ij' 


1 + Bi'Ki-ij + CiXi_ij'+i -I- AsXij'_i -I- RsX^j -I- CsXij'+i -I- AsXi+i^j-i -I- Bs'Ki+ij + CsXi+ij'+i + D 


0(B2) 
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where the 7 x 7-matrices Ai, Bi and Ci, for Z = 1, 2,3 are complex and defined in terms of the unperturbed stellar model. 

The elimination procedure that solves for the perturbation vectors Xij- is initiated by applying the boundary conditions 
at the stellar centre (* — 1 = 0), where we require and to vanish. The remaining five components of the perturbation 


vector Xo,j are located outside the grid, see Fig. (Al), and need not be defined. Hence, for i = 1 the terms with the coefficient 


matrices Ai, Bi and Ci disappear from the finite difference equation (FDE) (B2) for all j. 


By applying the boundary conditions at the stellar rotation axis {j — 1 = 0) and requiring and Fg to vanish, and noting 


that the remaining five components of the perturbation vectors Xi,o need not be defined (see Fig. Al), we conclude that for 


j — 1 the terms with the coefficient matrices Ai, with 1 — 1,2, 3, also disappear from the FDE (B2) for all i. 
For i = 1, j = 1 the simplified FDE’s can thus be written as 


(B3) 


B2X1U -|- C 2^1,2 + 733X2,1 + <73X2,2 + D — 0 

By inverting matrix B 2 we immediately find the relation: 

Xi,i = -B 2 "^C' 2 Xi ,2 - B^^B3^2,i - B2''C'3X2,2 - B^^D (B4) 

The elimination procedure is continued by stepping up j for given i = 1. For j = 2 the vector Xi,j_i no longer vanishes, 


but it can be eliminated by substituting relation (B4) for Xi,j_i into the FDE. After again multiplying by the inverted matrix 
BZ^ we obtain a relation of the form 


Xi,,' = a^’^Xi,,'+i + ^Pi’^X2,i -b 7 
1 = 1 


J + 1 




where a^’7, are complex 7 x 7-matrices and is a 7-vector. 


For j > 2 expression (B5) is applied after every step j to eliminate Xi,j from all previous steps (j — l,j — 2, 
repeating the elimination procedure until j = Ng we can apply the equatorial symmetry conditions 
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where <S is a diagonal 7 x 7-matrix. Its diagonal elements corresponding to and Fg are —1, while all remaining diagonal 
elements are 1. Using these symmetry relations to eliminate Xi,jvg+i, relation (B5) for Xi,iVj, can be re-expressed as 


Ng 


(B7) 


where and 7 ^’^® have been redefined. 

After substituting the above expression for Xi,jv^ in all previously obtained results, we have eliminated all unknowns 
Xi, k on the right hand side of expression (|B5[) and we end up with the set of relations (for k =1,2,3, .. Ng) 


Ng 


Xi,fc = ^/3/’''X2,,-by'’" 


(B8) 


where the and have been updated. 

We can now step up i and start a new strike j = 1, 2, ...Ng (for i = 2, 3, ...Ni). For i > 1 the perturbation vectors Xi_i,fc 
no longer vanish in the FDE. However, relations (B 8 ) obtained in the previous strike i — 1 can now be used to eliminate these 
perturbation vectors from the FDE. To this end the latter relations must be updated after every elimination step {i,j) to 
eliminate all Xi,; with I < j from the summation on the right hand side of (B 8 ). For step {i,j) this results in an expression 
of the form (for k = j — l,j,j + 1, ....Ng) 
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After substituting these results in the FDE (B2) for step {i,j), all terms left of 732 Xi,j -b ... can be eliminated and the FDE 
is cast in the form 
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where B 2 , C 2 and D have been redefined to collect the various coefficients of the corresponding vectors. 
By again multiplying the FDE with the inverted matrix the vector Xij' can be expressed as 


f+i 


i + l,l + 


Ng 
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l^j+2 




(B12) 


which is of the same format as the result of the previous elimination step (BIC). The current elimination step {i,j) is terminated 
by updating expressions (B9-B10). To this end expression (BB) for Xij is substituted on their right hand sides. Equation 
( [B9|) must be updated for all k > j. 

After repeating this elimination procedure until j = Ng the equatorial boundary conditions can be applied after which 
the stored relations for this strike i can be written as {k = 1, 2, ...Ng) 

Ng 

Xi.fe = ^/3;’''Xi+i,, + 7 '’'“ (B13) 


with which we can start a new strike i + 1 up to i = A^r, which corresponds to the stellar surface. At the surface we require 
the Lagrangian pressure perturbation 5P to vanish and the radial flux perturbation to fulfill Stefan-Boltzmann’s relation. In 
this way XjVr+i.j can be ignored, so that the first summation on the right hand side of relations (B9, BIO) and equation( |Bll| ) 
vanishes. After working through the last strike i = N^ we finally arrive at j = Ng. Note that now also the second summation 
in the above relations has vanished. Substituting the thus simplified relations (Bf-BIC) and applying the equatorial symmetry 
relations, all terms in the FDE can now be collected in B 2 and D only 


B2^Nr,Ng + D = 0 


(B14) 


Solving this equation we find the first perturbation vector X]v,,jve = —B^^D. We can find all other perturbation vectors by 
stepping back through the entire elimination procedure. E.g. having solved X]v,,jvs, X]v,,iVe-i, XjVr,j+i, we find Xw^.j by 
direct substitution into the stored relations (B12) 


X ... .j 

NrJ =0 


Ng 

XjV,,J + l + E 

l=j + 2 


Pi "’■'Xiv,,i +7 


(B15) 


Having solved in this way for all Xjv,,i, we find XjVr-i,i from relations (BIC), etc. 
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